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ABSTRACT 

Using three-dimensional non-equilibrium ionization (NEI) hydrodynamical simulation 
of the interstellar medium (ISM), we study the electron density, rig, in the Galactic 
disk and compare it with the values derived from dispersion measures towards pulsars 
with known distances located up to 200 pc on either side of the Galactic midplane. 

The simulation results, consistent with observations, can be summarized as follows: 
(i) the DMs in the simulated disk lie between the maximum and minimum observed 
values, (ii) the log(ne) derived from lines of sight crossing the simulated disk follows 
a Gaussian distribution centered at /i = —1.4 with a dispersion a = 0.21, thus, the 
Galactic midplane {Uf.) = 0.04 ± 0.01 cm^'^, (iii) the highest electron concentration 
by mass (up to 80%) is in the thermally unstable regime (200 < T < 10^-^ K), (iv) 
the volume occupation fraction of the warm ionized medium is 4.9-6%, and (v) the 
electrons have a clumpy distribution along the lines of sight. 
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1 INTRODUCTION 

Dispersion measures (DMs) towards pulsars with known dis- 
tance d can be used to derive the mean electron density in 
the Galaxy through the relation (ue) — DM/d. For pulsars 
located at |b| < 5° the derived mean electron density in the 
Galactic plane is {ue) ~ 0.02 - 0.1 cm"^, and 0.01 - 0.017 
cm~'^ in the spiral arms and inter-arm regions, respectively 
(Ferriere 2001; Gaensler et al. 2008). Berkhuijsen & Fletcher 
(2008), using the DMs towards 34 pulsars, mostly outside 
the galactic plane, showed that the log(ne) follows a Gaus- 
sian distribution, relating it to the nature of the interstellar 
turbulence (see e.g. reviews by Elmegreen & Scalo 2004). So 
far, numerical (magnetized and unmagnetized) 3D simula- 
tions have assumed colhsional ionization equilibrium (CIE) 
conditions for the ISM. In fact the thermal evolution of the 
ISM is determined by heating and cooling processes, which 
in general are not synchronized with ionization and recom- 
bination processes, respectively (e.g., Kafatos 1973; Shapiro 
& Moore 1976). Hence, below lO'' K, deviations from CIE 
conditions occur, thereby affecting the ionization structure 
of the interstellar gas, and thus the local electron density. In 
particular, if delayed recombination plays a role, the number 
of free electrons may be severely underestimated (Breitschw- 
erdt & Schmutzler 1994). 
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Here we study the electron density, rie, and it's mean, 
(rie), in the Galactic disk (up to 200 pc on either side of 
the midplane) using the first to date global hydrodynamical 
simulation of the interstellar gas, evolving under NEI con- 
ditions. Furthermore, we compare simulation results with 
estimates of (rie) obtained from available pulsar DMs. In 
forthcoming papers we explore the rie distribution and it's 
topology in the thick disk and halo of the Milky Way and 
other galaxies. 

This paper is organized as follows: Section 2 deals with 
the model setup; Sections 3 and 4 present the observational 
data, and simulation results, respectively. A discussion and 
final remarks in Section 5 close the paper. 



2 MODEL AND NUMERICAL SETUP 

We simulate hydrodynamically the supernova-driven ISM in 
a patch of the Galaxy centered at the Solar circle with an 
area of 1 kpc^ and extending to 15 kpc on either side of 
the Galactic midplane following de Avillez & Breitschwerdt 
(2005, 2007; AB0507). The simulations are carried out with 
the EAF-parallel adaptive mesh refinement code coupled to 
the newly developed E(A+M)PEC cod«Q Avillez & Bre- 

1 www.lca.uevora.pt/research.html presents a description of the 
code, ionization fractions, cooling, and emission spectra. 
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itschwerdt 2012) featuring the time-dependent calculation 
on the spot (at each grid cell) of the ionization structure of 
H, He, C, N, O, Ne, Mg, Si, S and Fe and emissivities. 

The physical model includes supernovae (SNe) types la, 
Ib-f c, and II, a gravitational field provided by the stars in 
the disk, local self-gravity (excluding the contribution from 
the newly formed stars), heat conduction (Dalton & Balbus 
1993), uniform heating due to a UV radiation field normal- 
ized to the Galactic value and varying with z, and photoelec- 
tric heating of grains and polycyclic aromatic hydrocarbons. 

E(A-|-M)PEC uses the recommended abundances of As- 
plund et al. (2009), and calculates electron impact ioniza- 
tion, inner-shell excitation auto-ionization, radiative and di- 
electronic recombination, charge-exchange reactions (recom- 
bination with Hi and Hel , and ionization with Hll and Hen ), 
continuum (bremsstrahlung, free-bound, two-photon) and 
line (permitted, semi-forbidden and forbidden) emission in 
the range lA -610/i. The code also includes ionization of 
Hi by Lyman continuum photons emitted during the recom- 
bination of helium. The internal energy of the plasma in- 
cludes the contributions due to the thermal translational 
energy plus the energy stored in (or delivered) from high 
ionization stages. Electron impact ionization rates are taken 
from Dere (2007), while radiative and dielectronic recom- 
bination rates are based on AUTOSTRUCTURE calcula- 
tions (Badnell et al. 2003, BadneU 200680 including the 
latest corrections to Fe ions by Badnell (2006b), Nikolic et 
al. (2010), and Schmidt et al. (2008). Radiative and dielec- 
tronic recombination rates for Slland Fevilare from Maz- 
zotta el al. (1998). For the remaining ions we adopt the total 
recombination rates derived with the unified electron-ion re- 
combination method (Nahar & Pradhan 1994) and available 
at NORAD-Atomic-Dat4!l 

A coarse grid resolution of 8 pc is used, while the finest 
AMR resolution is 0.5 pc (4 levels of refinement) for \z\ 5C 2 
kpc, 4 pc for |z| > 4 kpc, and 1 pc elsewhere. Periodic and 
outfiow boundary conditions are set along the vertical faces 
and top/bottom {z = ±15 kpc) of the grid, respectively. 



3 THE ELECTRON DENSITY IN THE DISK 

To make a detailed comparison between simulation results 
(discussed below) and the electron density distribution de- 
rived from DMs towards pulsars, we reassessed the existing 
data to select pulsars with best possible distance estimates. 
Density PDFs were created using pulsars located up to 8 kpc 
away from the Sun and with \z\ < 200 pc from the Galaxy's 
midplane. 

For our analysis we chose pulsars from the ATNF cata- 
with independent distance esti- 
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mates, i.e., estimates without using NE2001 model (Cordes 
& Lazio 2003; CL03). These estimates resulted from parallax 
measurements, absorption-line (Hl21-cm or OH line) stud- 
ies, physical associations or timing towards pulsars in the 
galactic disk {\z\ < 200 pc). Most parallax estimates were 
obtained from Verbiest et al. (20101, which are updated by 
Chatterjee, SQ 



^ amdpp.phys.strath.ac.uk/tamoc/DATA/ 

3 www.astronomy.ohio-state.edu/^nahar and references therein. 
^ [iiJvww.astro.cornell.edu/ shami/psrvlb/parallax.html" 
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Figure 1. Histogram (solid black line) and Gaussian fit (dashed 
black line) of log(?ie) distribution as obtained with a series of 
24 pulsars located at 200 < distance < 8000 pc from the Sun, 
with \z\ ^ 200 pc. The mean value is Iog(ne> = -1.47 ±0.02 and 
cr = 0.17 ± 0.02. Red solid and dashed lines refer, respectively, to 
the histogram and Gaussian fit (centered at fi = —1.44 cm""^ and 
having a dispersion a = 0.21) of the ('T'e}35o_4oo Myr (§^)- 



For our study we augmented the ATNF catalog by in- 
cluding distance estimates from absorption studies wherever 
available from the literature, including two measurements 
from associations (Frail et al, 1996; Saravanan et al. 1996; 
Johnston et al. 2003; Minter et al 2008; Weisberg et al 2008; 
see references therein). To quantify errors in distance esti- 
mates, we compute the fractional error involved in parallax 
and absorption studies from the difference between the up- 
per and lower estimates. Firm error estimates are not avail- 
able for several pulsars in the ATNF catalog, whereby the 
fractional error was computed from the difference between 
the distance estimate and the best estimate available from 
NE2001 model (parameter DISTl in ATNF catalog). Errors 
in distance estimates using parallax measurements range be- 
tween 2 — 35%. The median error in distance estimate in our 
overall sample is 18%, a significant improvement over esti- 
mates using NE2001 model (CL03). 

A total of 122 pulsars with reliable distance estimates 
were thus selected, 33 of which have \z\ < 200 pc. Nine 
pulsars were removed from this sample as they are located 
within the vicinity of the Sun, and measurements of elec- 
tron density towards those sightlines are significantly af- 
fected by local structures, such as the Local Bubble, Gum 
nebula and Loop I (e.g., CL03). The histogram of the com- 
puted log{ne) derived from the sample DMs and the best-fit 
gaussian curve are shown in Figure [l] (black lines). The fit 
is centered at log{ne) = —1.47 ± 0.02 and has a dispersion 
of cr = 0.17 ±0.02, with a x^/dof = 1.62. These values com- 
pare well with those derived directly from our pulsar sample 
data: log(ne) ~ 1.42 and a ~ 0.23. 



4 SIMULATION RESULTS 

The simulated supernova-driven interstellar medium is char- 
acterized by several evolutionary phases that have already 
been described in previous works: (i) domination of initial 
conditions being only wiped out after some 80 Myr of evolu- 
tion, (ii) the full establishment of the continuous disk-halo- 
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Figure 2. Total (left) and electron (right) density distributions (in log scale) in the Galactic midplane at 400 Myr of evolution. Red 
regions in the left panel represent high density material, with molecular clouds being represented by black. Electron densities smaller 
than log rig = —4 (white regions; right panel) are located in both high (atomic and molecular clouds) and lower (bubbles) density regions. 



disk circulation (also known as the Galactic fountain), and 
(ill) the dynamical equilibrium (occurring at 200-220 Myr 
of evolution) in a statistical sense that determines the dy- 
namics of the interstellar medium in the Galactic disk and 
its interaction with the halo. As a result of this evolutionary 
path driven mainly by supervovae, the ISM becomes frothy 
and turbulent with a mean Mach number of 3 (AB0507). 
Low temperature gas (high density) is concentrated into fil- 
amentary structures and molecular clouds (black regions in 
the left panel of Figure |2] showing the Galactic midplane 
density at 400 Myr of evolution), while hot gas (low density 
regions in the same panel) is concentrated into bubbles and 
superbubbles that dominate the landscape. 

The turbulent nature of the simulated ISM, the on- 
going physical processes and driving mechanisms (e.g., SNe 
and stellar winds) affect the electron distribution, which fol- 
lows the topology of the medium. Most of the electrons are 
distributed into the thermally unstable regime (200 < T < 
-|^q3.9 j^-j ^ volume (/„) and mass (/a/) filling fractions 
of 56 — 60% and 77—80%, respectively (Figure |3] top panel); 
32 — 40% of the electron mass is locked in filamentary struc- 
tures and shells, having /„ = 6 — 7%. The electron distri- 
bution in warm ionized medium has fv — 4.9 — 6% and 
/m = 7 — 10%. These are time-dependent variations whose 
time average (over a period of 100 Myr using 1001 disk snap- 
shots taken at every 0.1 Myr) histogram of log He (Figure |3] 
bottom panel), is fitted with a composition of two gaus- 
sians centered at logrie = —1.437 (for T< 10^'^ K gas) and 
logrie = —2.64 (for T> 10^'^ K gas), respectively, with the 
latter representing mostly photoionized regions. 

Further insight into the simulated electron distribution 
can be obtained directly from DMs along LOS, with length 
d, crossing the simulated ISM, allowing the determination of 
the mean electron density. The LOS have increasing lengths 
(with a step length of 10 pc) up to 1 kpc from a vantage 
point located at {x, y, z) — (0, 0, 0) (left bottom corner of 
the electron density maps in Figure |2| for l^j ^ 200 pc and 
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Figure 3. Top panel: Time evolution between 300 and 400 Myr 
of the electrons volume (solid lines) and mass (dashed lines) oc- 
cupation fractions in the simulated disk for thermally unstable 
regime (orange), gas in shells and filamentary structures (red), 
and warm ionized medium (WIM; cyan). Bottom panel: Time 
averaged histogram of the total (black ), and the two gaussian 
fits (dashed black lines) of the electron density in the simulated 
disk over the evolution time 300-400 Myr. The electron density 
in different temperature regimes are displayed. 
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Figure 4. Angular distribution of tlie DMs along LOS witii 
lengths varying between 100 and 1000 pc at 350 and 400 Myr 
(top two panels) and tiie time averaged {DM)t (bottom panel) 
over the evolution time 350-400 Myr). 



spanning 90° (with a 1° separation) when projected onto 
the Galactic midplane. 

Figure [4] shows the DM along the LOS with lengths 
varying between 100 and 1000 pc, with a step length of 100 
pc at 350 and 400 Myr (top two panels) and its time average 
(bottom panel) over 501 disk snapshots taken every 0.1 Myr 
for 350 — 400 Myr of evolution. Due to the ongoing turbulent 
processes in the Galactic disk the DM has variability along 
and between the lines of sight (Figure[4] top two panels). The 
longitudinal variation is due to the inhomogeineity of the 
ISM (e.g., bubbles and superbubbles) as well as due to the 
turbulent nature of the regions crossed by the LOS. These 
variations with distance and among LOS are still present in 
the time averaged DM, although the longitudinal variations 
are smoothed in the averaging process. However, a record of 
any event that prevailed over a large period of time (e.g., 
a superbubble whose contribution to the DM is small as in 
the case of the observed hump in the DM between 48° and 
65°) is kept. 

We further compare the DMs and (rie) derived from 
pulsars observations (§3), with those resulting from averag- 
ing over all the LOS crossing the simulated disk {\z\ < 200 
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Figure 5. Minimum (dotted line), mean (solid line) and maxi- 
mum (dashed line) of the averaged DM (top panel) and rie (bot- 
tom panel) for all LOS crossing the disk up to 4 kpc from the 
Sun at 310 Myr ({X)^-^^^ ^^j.; red lines), 360 Myr {(X)^^^ j^yj.; 
green lines), and over the period 350-400 Myr ((^)35o_4oo Myr' 
black lines). The blue circles represent the DMs (top) and (we) 
(bottom) derived from pulsars observations. 



pc) and up to 4 kpc from the Sun at specific times (e.g., 310 
and 360 Myr) and over the period 350-400 Myr. Top panel 
of Figure [5] displays the observationally derived DMs and 
the averaged DMs at 310 Myr {{DM)^-^^ Myr' ''^'^ hues), 
360 Myr {{DM)^^^ Myr' grssn lines), and over the period 
350-400 Myr {{DM)^^g_^gf^ jyjyj.; black hues). 

Due to the limitation of the box size in our simula- 
tion, we calculated the DMs for distance d > 1 kpc taking 
advantage of the periodic boundary conditions used in the 
simulations. The DMs were calculated from different van- 
tage points in the midplane (including (x = 0,y = 0,x = 0) 
and between 500 and 1400 pc from it) using LOS crossing 
the disk volume at different angles, up to 4 kpc from the 
vantage point and up to z = ±200 pc. This procedure, al- 
though inevitable here, may introduce artefacts in the sim- 
ulation results, but we do not expect that our results are 
strongly affected though, because the simulation shows a 
certain pattern repetition in the ISM density and tempera- 
ture distributions on scales of the order of a few correlation 
lengths (that is, a few times 75 pc, according to AB0507). 

The averaged (_DM)gg(j .^pj, lyfyj., after a steep growth 
in the first 100 pc, reaches a smooth increase for d > 500 pc, 
with the deviations of the minimum and maximum regard- 
ing the mean becoming constant over distance - consistent 
with the DMs variations in d observed in Figure [4] The aver- 
aged DMs at 310 and 360 Myr vary widely, but with values 
within the pulsar derived DMs. Therefore, the terms mini- 
mum and maximum of the (DM) Myr should not be 
taken as strict upper and lower values for all times, as they 
only represent an average over a specific time window of the 
simulation. 

After the large variation in the first few hundred par- 
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sees, the mean of the simulated {n^), for the cases shown 
in the bottom panel of Figure [S] reach a constant value 
(0.04±0.002 cm"^) and dispersion over large distances. This 
reflects the clumpy nature of the electrons in the turbulent 
ISM as would be expected if they are predominantly found 
in filaments and shells (Berkhuijsen & Miiller 2008). The 
histogram of the \og{ne) .^^^^^^g^ y[yj. and it's best Gaussian 
fit are displayed in Figure [T] by solid and dashed red lines, 
respectively. The fit is centered at log{ne) = — 1.4±0.01 with 
a = 0.2f ± 0.01. In the simulated disk, the clumpy material 
dominates the LOS and volume averaged PDFs which ex- 
plains the similar means between the electron density of the 
T< lO'*'^ K regime (Figure [sjl and the time averaged {ne)t 
(Figure [l] red lines). These fit parameters are similar to 
those displayed in Figure [l] (black lines) for observationally 
derived (rie). This similarity is indicative of the smoothing 
out of low and high density regions along the lines of sight, 
thus, not appearing in the PDFs, and leading to reduced 
Gaussian dispersion. 



5 DISCUSSION AND FINAL REMARKS 

In this letter we discuss the electron density distribution 
in the Galactic disk using the first to date 3D high resolu- 
tion NEI simulations of the ISM. The simulations trace the 
dynamical and thermal evolution of the interstellar gas, cal- 
culating on-the-spot the ionization structure, electron dis- 
tribution and cooling function of the gas at each cell of 
the grid into which the computational domain is discretized. 
We did not take into account stellar ionizing photons and 
their transport, but an averaged diffuse photon field, which 
should be a reasonable approximation, when describing the 
mesoscale ISM. The NEI structure modifies the cooling func- 
tion, which in turn enhances the number of free electrons. 

Both the simulated DMs and electron density are con- 
sistent with the observations: (i)nmost DMs lie within the 
maximum and minimum observed values, (ii) log{ne) is con- 
sistent with a gaussian distribution, (iii) the mean electron 
density is 0.04 ±0.01 cm~'^, and (iv) the volume filling frac- 
tion of the WIM is bracketed between 4.9 and 6%. These 
results lye in between the estimates by Gaensler et al. (2008) 
and Berkuijsen & Miiller (2008) for the electron density in 
the disk and WIM. Furthermore, the observed Gaussian 
distribution is present, irrespective of whether the ISM is 
isothermal (see discussion in Berkuijsen & Fletcher (2008) 
and refereces therein) or not. This can be explained in two 
complementary ways: a fully developed turbulent system 
can be considered as a large number of independent random 
variables, with the logarithm of each having a certain distri- 
bution approaching a Gaussian as the number of variables 
goes to infinity (Central Limit Theorem), being the case 
of supernova-driven ISM simulations reaching a statistical 
equilibrium (see, e.g., Vazquez-Semadeni & Garcia 2001). 
Or, alternatively, we can use the principle of maximum en- 
tropy. The information entropy of a continuous random vari- 
able X with probability density function p{x) is defined as 
H{X) = J ^^p{x) log p{x)dx. The lognormal distribution of 
X maximizes H{X), implying the least prior knowledge of 
the system (Sveshnikov, 1968). This is exactly what is ex- 
pected in homogeneous and isotropic turbulence with a large 



number of independent random variables. We study devia- 
tions from lognormal distributions in a forthcoming paper. 
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